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I/"") . Abstract 



We propose a new monotonically convergent algorithm which can 
enforce spectral constraints on the control field (and extends to arbi- 
trary filters). The procedure differs from standard algorithms in that 
at each iteration the control field is taken as a linear combination of 
■ the control field (computed by the standard algorithm) and the filtered 

£ ' field. The parameter of the linear combination is chosen to respect the 

monotonic behavior of the algorithm and to be as close to the filtered 
field as possible. We test the efficiency of this method on molecular 
alignment. Using band-pass filters, we show how to select particu- 
lar rotational transitions to reach high alignment efficiency. We also 
consider spectral constraints corresponding to experimental conditions 
using pulse shaping techniques. We determine an optimal solution that 
could be implemented experimentally with this technique. 



■ 1 Introduction 



Quantum control is a field of growing interest both from the experimental 
and theoretical points of views [D [3]. On the experimental side, the 
development of pulse shaping techniques opens new ways to control atomic 
or molecular processes by laser fields. Many promising results have been 
obtained in the implementation of closed-loop control experiments (CLC) 
based on genetic or evolutionary algorithms (See [H El EJ [TJ [8] to cite a 
few). This setup is made of a pulse shaper controlled by the algorithm 
which from the results of the preceding experiments builds an improved 
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new control field. Such algorithms lead to very efficient solutions, but have 
some negative points. In particular, no insight into the control mechanism 
is gained from this approach since no knowledge (or a minimal one) about 
the system is needed and the control field is not optimal by construction. 

On the theoretical side, optimal control theory (OCT) is a powerful tool 
to design electric fields to control quantum dynamics [9j [TUl fTT] . Numerical 
aspects of optimal control theory have been largely explored. For simple 
systems with few energy levels, geometric aspects of optimal control based 
on the Pontryagin maximum principle \1'2\ 113] have also been investigated 
[14|, \TE[ 116] . Monotonically convergent algorithms is other efficient approach 
to solve the optimality equations [13 [HI [El EQl EU [22l [231 E31 [^1 [26]. 
They have been applied with success to a large number of controlled quan- 
tum systems in atomic or molecular physics and in quantum computing. 
These methods are flexible and can be adapted to different non standard 
situations encountered in the control of molecular processes. Among recent 
developments, we can cite the question of nonlinear interaction with the con- 
trol field [271 ESI [291 EQ] and the question of spectral constraints on the field 
[3"T| [3"2"1 1551 EH E2] . This latter problem is particularly important in view 
of experimental applications since not every control field can be produced 
by pulse shaping techniques [H El El E] • For instance, liquid crystal pulse 
shapers are able to tailor only a piecewise constant Fourier transform of the 
control field in phase and in amplitude. Experimentally, the spectral am- 
plitude and phase are discretized, e.g. into 640 points which is the number 
of pixels in a currently used standard mask. In this context, a challenging 
question is to incorporate advantages and efficiency of OCT into CLC. OCT 
should provide insights into the control mechanisms and speed up the con- 
vergence of the experimental algorithm towards an optimal solution. Note 
that the theoretical analysis of the control dynamics is more efficient if the 
computational schemes can include experimental constraints directly in the 
algorithm. A brute force strategy which consists in applying the constraints 
after the optimization leads to worse results. 

This paper aims at taking a step towards the answer to this question. 
For that purpose, we present a monotonically convergent algorithm which 
can take into account spectral constraints on the control field. Unlike other 
proposals in the literature [311 [32J, [331 EH ES] our approach exactly enforces 
the monotonicity property which is important to ensure the convergence of 
the algorithm. The construction of a monotonic algorithm with spectral 
constraints is a very difficult task, since constraints in the frequency domain 
require nonlocal knowledge of the control field at all times whereas this field 
is computed progressively in time by the algorithm. These two requirements 
are incompatible. This difficulty has been bypassed in |31| for finite impulse 
response filters by using the control field of the preceding iteration to con- 
struct the Fourier transform in order to satisfy spectral constraints. Here we 
propose a similar algorithm but with a modified electric field. This means 
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that we consider a standard monotonically convergent algorithm [19} [20| [2T] 
giving at each iteration a control field. In particular, only the energy of 
the field is penalized by the cost which does not depend on spectral con- 
straints. We then apply a filter to this field to obtain a filtered control field. 
We next construct a new field as a linear combination of these two fields. 
The parameter of the linear combination is such that the algorithm remains 
strictly monotonic and the new control field is as close to the filtered field 
as possible. As expected, the filtering comes with the price of a slow down 
in convergence. 

We test the efficiency of this approach on molecular alignment which is a 
well-established topic in quantum control [361 EH EE] ■ Molecular alignment 
has been optimized both experimentally [391 HOI Wl\ [4"2] and theoretically 
in different studies using genetic algorithms |43[ [4"4"1 [4"5] or optimal control 
theory [461 HT] . Note that previous optimal control studies have not consid- 
ered spectral constraints. The use of spectral constraints will allow one to 
construct new optimal control fields to reach a high alignment efficiency In 
[46j , it was shown that aligned states are reached by rotational ladder climb- 
ing, i.e., by successive rotational excitations. Using spectral constraints, we 
determine an optimal control field that induces only particular rotational 
transitions. This shows how to guide the algorithm towards a particular 
mechanism. In a second step, we consider spectral constraints correspond- 
ing to experimental pulse shaping techniques. Our new algorithm leads to 
an optimal solution that could be implemented experimentally. The solu- 
tions determined from genetic algorithms in [43 \ [4"4"1 14"5] do not correspond 
to an optimal solution both maximizing the alignment and minimizing the 
energy of the control field. Minimizing the intensity of the field can be in- 
teresting to avoid parasitic phenomena such as ionization. Finally, we recall 
that in molecular alignment, due to the rapid oscillations of the laser field, 
the effect of the permanent dipole moment averages to zero and plays no role 
in the control of the dynamics. The interaction between the molecule and 
the electric field is therefore of order 2. Since the interaction is nonlinear, 
we use the monotonic algorithm introduced in [27] with a non-standard cost 
of power 4 in the electric field to compute the optimal solution. 

The paper is organized as follows. In Sec. [21 we present the problem of 
laser induced alignment of a linear molecule by non-resonant laser fields. We 
describe how monotonically convergent algorithms can be applied to such 
systems. We explain how to modify the algorithm to take into account spec- 
tral constraints. We show the efficiency of this new approach in Sec. [3] for a 
filtering which corresponds either to three particular rotational frequencies 
or to experimental conditions using pulse shaping techniques. Some conclu- 
sions and discussions are presented in Sec. [U Some technical computations 
are reported in appendix |A"1 
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2 Optimal control of molecular alignment with spec- 
tral constraints 



2.1 Description of the model 

We consider the control of a linear molecule by a non-resonant linearly po- 
larized laser field of the form £(t) = E(t) cos(u;t) where oj is the carrier wave 
frequency and E(t) the laser pulse envelope [361 EH EH]- In the case of a 
zero rotational temperature, the dynamics of the system is governed by the 
Schrodinger equation. In a high-frequency approximation [36] , this equation 
can be written as 

(1) *J^(*)> = [BJ 2 -\E(t) 2 (Aacos 2 + a± )]\m) 

where B is the rotational constant and Aa = an — a± is the difference be- 
tween the parallel and perpendicular components of the polarizability tensor. 
We use atomic units unless otherwise specified. For numerical applications 
we consider the molecule CO with the following parameters B = 1.931 cm -1 , 
an = 15.65 and qj_ = 11.73 in atomic units [48] . 

2.2 Monotonically convergent algorithm 

The optimal control problem is solved by a monotonically convergent algo- 
rithm. The goal of the control is to maximize the projection of the system at 
time tf onto a target state where tt is the duration of the control. The tar- 
get state is taken as the state which maximizes the alignment in a reduced 
finite-dimensional Hilbert space spanned by the j op t first rotational levels 
[Ml [32 EI] . The reduced Hilbert space is denoted ~Hj opt - We assume in this 
paper that j op t = 8. A basis of the Hilbert space is given by the spherical 
harmonics \j, m) with j > and —j<m<j. In case of pure state systems 
at a temperature T = 0K, the target state \cj)f) is simply the eigenvector of 
maximum eigenvalue of the projection of the operator cos 2 9 onto T~Lj opt ■ The 
optimal control of molecular orientation and alignment has been considered 
in a series of articles [Ml E7J EH [27] . The novelty of our approach lies in 
the fact that we will consider in addition spectral constraints which leads 
to new optimal control fields. In other words, this means that there exists 
no unique optimal control field to reach a given target state and that it is 
possible to select an optimal solution via spectral constraints. Due to the 
interaction of order 2 between the system and the field, we use the algorithm 
introduced in [27] with a cost which is quartic in the field. 

We consider the initial state |0o) = |0, 0) and the following cost func- 
tional: 

(2) J=\^ f m f ))\ 2 - f f X{t)E{tfdt 

Jo 
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where X(t) = Ao/s(t) is a positive function. Following [53], we will choose 
s(t) = s'm 2 (7rt/tf) which penalizes the amplitude of the pulse at the be- 
ginning and at the end of the control. Here Ao is chosen to express the 
relative weight between the projection and the energy of the electric field. 
The augmented cost functional J is defined as follows 

J = |»(t / ))| 2 -^ / A(i)^(t) 4 dt-2Im[(^(i / )|0 / ) j\ X it)\^-H)\m)dt] 

where \x(t)) is the adjoint state and Im denotes the imaginary part. Setting 
the variations of J with respect to \ip(t)}, \x(t)) and E(t) to be zero implies 
that \ip(t)) and \x(t)) satisfy the Schrodinger equation with respective initial 
conditions |</>o) and \4>f)'- 

(»§-#(*)) hKt)> = o 

(4-H(t))\x(t))=0 ■ 

The optimal electric field is solution of the polynomial equation: 

(4)4A(t)S(t) 3 +Im[(^ / )| < /. / }(x(i)l(Aacos 2 + Q ± ) J E(t)|^(t))] =0. 

We solve this set of coupled equations by using a monotonic algorithm. 
To simplify the presentation of the algorithm [27] . we consider that the 
forward and backward electric fields are the same. The iteration is initiated 
by a trial field Eo(t). Let us assume that at step k of the iterative algorithm 
the system is described by the triplet (\ipk(t)}, |Xfc-i(*))> E^it)) associated 
to the cost Jfc given by 



(5) Jk = J(E k ) = \(<l>f\il>k(tf))\- Kt)E k (t)dt. 

Jo 

The spectral constraints are defined through a filter T in the frequency- 
domain. This means that for a given electric field E(t), J-[E(t)] only contains 
the admissible frequency components of the field. 

We determine the triplet at step (k + 1) from the triplet at step k by the 
following operations. We first propagate backward the adjoint state \xk(t)) 
with the field Ek{t) and initial condition | </>/): 

(6) i^\ Xk (t))=H(E k (t))\x k (t)). 

We then propagate forward the state \tpk+i) with initial condition |0o): 

(7) ij^fc+i(*)> = H(E k+1 (t))\^k+i(t)), 

computing at the same time the electric field Ek+\(t) by the procedure 
explained in the appendix [XI 
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Having computed the non-filtered optimal field E k +i(t), the next step of 
the algorithm consists in introducing the filtered field E k+ \ as: 

( 8 ) E k+ i{t) = E k+ i^ k+1 (t) 
where we have denoted for any /i G [0, 1]: 

(9) E k +iAt) = » E k+i(t) + (1 - n)F[Ew(t)] 

and Hk+i is such that the algorithm remains monotonic, i.e. such that 
AJ k+1 (fj, k+1 ) = J(E k+1;IJ , k+1 ) - J{E k ) > 0. We choose a value of fj, k+1 such 
that A Jfe+i(/^fc+i) > but close to 0. We explain in Sec. [3] on the example 
of molecular alignment how to determine numerically fik+i- This choice 
ensures the monotonic behavior of the algorithm and to be as close to the 
filtered field as possible. 

3 Numerical results 

We first present numerical results at T = K to show the efficiency of the 
algorithm for long control fields (Sec. 13. ip . Spectral constraints correspond 
in this case to band-pass filters which allow to select particular rotational 
transitions. We also analyze the temperature effects for shorter control 
fields in Sec. 13.21 with constraints mimicking the experimental pulse shaping 
techniques. 

3.1 Zero temperature 

We consider a control field with a duration of ten rotational periods. The CO 
rotational period is denoted by t per . This long duration allows to simplify 
the structure of the Fourier transform of the fields. The filter is composed 
of three bandpass filters around the frequencies AB, 10B and 26B with a 
bandwidth of B/2. Combining these three frequency components, we obtain 
the different rotational transitions between the even rotational levels which 
are populated during the control. We have for instance that W02 = WB — 4B, 
W24 = 10JB + 45 and uj^q = 26B — 4B where ujjji is the difference of 
energy between the Jth and the J'th rotational levels [29J. This shows 
that all the even rotational levels can be populated using only these three 
frequency components. Also, this means that the target state \4>f) can be 
reached from the initial state \4>o)- To compute the parameter fj, k +i at each 
iteration, we have used a dichotomy method to determine the zero, fj,Q, of 
AJfc +1 (/i) = J(Ek+i lft ) — J(E k ). The value fJ> k +i is the first value of /i given 
by this approach such that AJfc + i(/x) > and \fi — fio\ < 0.01. As could be 
expected, the smaller the difference \fi — fj>o\ is, the slower the convergence 
of the algorithm is. 
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We apply the standard algorithm |27J and the one with filtering to max- 
imize the projection onto the target state at time t = tt. The trial field is 
a Gaussian pulse with a duration (FWHM) equal to 3 ps. The parameters 
of the algorithms are taken to be Ao = 1 and rj = 1. The results of the 
computations, i.e. the evolution of (cos 2 0)(t) and the optimal control field, 
are presented in Fig. [TJ for the algorithm with filtering. Very good results 
are obtained with a final projection \{ift(tf)\4>f)\ 2 larger than 0.99. Similar 
efficiencies have been reached by the standard optimal field displayed in Fig. 
[Tb. As could be expected, the filtering slows down the convergence of the 
algorithm. We get a projection close to 0.98 after 200 and 400 iterations 
for respectively the standard and the modified algorithms. The monotonic 
behavior of the new algorithm can be seen in Fig. UJl. The evolution of the 
parameter [i^ as a function of the number of iterations is also represented 
in Fig. [JJi. One sees that is equal to zero for the first 300 iterations and 
then it increases to reach a value close to 1 for k — 600. Since after 300 
iterations, the optimal field for \i = 1 and the completely filtered one for 
fi = are very close to each other, the non-zero value of [i^ does not affect 
the filtering of the algorithm. The discontinuous behavior of /i^ in Fig. [TH 
is due to the dichotomy method and to the criterion used to determine this 
parameter. Note that the evolution of is the same in all the examples 
we have considered, i.e. a zero value for the first iterations and then a slow 
increase towards the value 1. Figures [2] show the Fourier transform of the 
optimal solutions obtained by the standard algorithm and the one with fil- 
tering. The spectral structure of the filtered optimal solution is very simple 
with respect to the one without filtering. As expected, only three frequency 
components appear in the Fourier transform of the filtered field, whereas no 
frequency component can be clearly distinguished for the optimal solution. 
One sees that the increase of the parameter ji^ after 400 iterations does not 
affect the Fourier structure of the optimal solution. Figure [Tb shows that 
the standard optimal solution does not use the ten rotational periods since 
the control field vanishes after 1.75 periods. This short duration leads to a 
very complicated spectral structure in Fig. [2b which contrasts with the very 
simple structure of the filtered optimal field. The spectral constraints can 
be viewed here as a tool for guiding the control dynamics along a selected 
pathway in the frequency domain in the sense that the standard optimal so- 
lution is a superposition of many pathways leading to a complicated Fourier 
spectrum. 

3.2 Non-zero temperature 

In this section, we test the efficiency of the algorithm at a non-zero tem- 
perature. The system is described by a density matrix with a dynamics 
governed by the von Neumann equation [371 [38] . The initial density matrix 
is the Boltzmann equilibrium density operator at temperature T. We use 
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Figure 1: Plot as a function of the adimensional time t/t per of (a) the 
expectation value (cos 2 6) obtained from the algorithm with filtering, (b) 
the corresponding optimal filtered field and (c) the standard optimal field, 
(d) Plot as a function of the number of iterations k of the adimensional cost 
Jk (solid line) and of the parameter (dashed line) for the algorithm with 
filtering. 
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Figure 2: Normalized square modulus of the Fourier transform of the trial 
field (a), of the optimal field obtained by the standard algorithm (b) and 
of the optimal field obtained by the algorithm with filtering. The vertical 
dashed lines correspond to the rotational frequencies between two successive 
rotational levels. 
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the target state p pt constructed in Ref. [50] for a linearly polarized elec- 
tric field. This density matrix acts in the space 7ij t with j op t = 8 which is 
chosen with respect to the temperature and the intensity of the field used. 
As before, we have used the algorithm of [27\ with a cost functional similar 
to the one of Sec. [2J We have considered three temperatures T = 5, 7 
and 10 K and a shorter control field with a duration equal to one rotational 
period. The numerical parameters of the trial gaussian field are taken to be 
37.5 TW/cm 2 for the peak intensity and 1 ps for the duration (FWHM), 
except for T = 7 K where the duration is equal to 0.475 ps. 6000 iterations 
have been used for the first two temperatures and 8000 for the last one. The 
parameters of the algorithm are equal to Ao = 1 and rj = 1. 

We have used a filtering which mimics experimental control techniques 
using liquid crystal pulse shapers |44J. Such pulse shapers work in the fre- 
quency domain by tailoring the spectral phase and amplitude of the electric 
field. This operation can be defined through a filter T which transforms over 
a given bandwidth the Fourier transform of the control field into a piecewise 
constant function in phase and in amplitude. The bandwidth is taken to 
be 7.28 THz which corresponds to two times the 5th rotational frequency. 
A larger bandwidth has given equivalent results. The filtered control field 
J-{E(t)) is finally obtained by an inverse Fourier transform. In the following 
computations, the number of pixels N is equal to 64, 128 or 256 both in 
amplitude and in phase. The pixels that discretize the Fourier transform 
are taken equally spaced in the frequency interval defined by the bandwidth. 

To determine the parameter fj,k+i, we compute at each iteration AJ/t +1 (/i) = 
J(Ek+i^) — Jk(Ek) for 10 values of \i equally spaced in the interval [0, 1]. 
Using these ten values, we construct a polynomial fit of A J/ s+ i( / u). We 
define Hk+i as the smallest value of \i such that AJfc +1 (yUfc +1 ) = 0.01 x 
max /ie[o,i] [A Jfc+i (/x)]. We have checked that these ten values are sufficient. 

The results obtained for the three temperatures are displayed in Fig. 
with a fixed number of pixels equal to 128. Very efficient optimal fields have 
been constructed by the algorithm since a final projection of 99.2%, 97% 
and 93.5% has been obtained for respectively T = 5, 7 and 10 K. As can 
be expected, the temperature has a negative effect on the alignment but 
the control remains robust with respect to the temperature. Note also the 
completely different structures of the three optimal fields leading to three 
different evolutions for the projection. This example shows that our algo- 
rithm can be used to construct realistic control fields both in temperature 
and in intensity. In Fig. [3b, one observes a similar evolution for the three 
parameters jik which are close to zero for the first hundred iterations and 
then quickly increase for k > 100. The increase of /i& is more pronounced 
for higher temperatures which indicates that spectral constraints are more 
difficult to impose when the temperature increases. 

Figures [U illustrate the impact on the algorithm of the number of pixels. 
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Figure 3: Plot as a function of the adimensional time t/tf of (a) the optimal 
electric field and (b) the projection onto the target state p op t- (c) Plot as 
a function of the number of iteration k of the parameter Solid, dashed 
and dotted lines correspond respectively to T = 10, 7 and 5 K. The number 
of pixels is equal to 128. 






600 1200 1800 2400 3000 

k 



Figure 4: Same as Fig. [3] but for T = 5 K and different numbers of pixels. 
Solid, dashed and dotted lines correspond respectively to N = 128, 64 and 
256 pixels. 



11 



We have obtained a final projection of respectively 0.975, 0.989 and 0.992 for 
N = 64, 128 and 256 pixels. In Fig. [H we have represented the fields after a 
final filtering to get solutions that can be implemented experimentally. We 
then obtain 0.813, 0.892 and 0.926 for the three cases. The final result passes 
for N = 64 from 0.975 to 0.813 which shows that the optimal solution does 
not satisfy very well the spectral constraints. One deduces that N = 64 
is not a sufficient number of pixels to construct the optimal field. This 
observation is also related to the rapid increase of /i£ in Fig. Eh for N = 64. 
To summarize, this result indicates the minimum number of pixels that has 
to be used experimentally to control molecular alignment with an optimal 
control field. 

4 Conclusion 

We have proposed a new monotonically convergent algorithm to take into 
account spectral constraints. The procedure is built on the standard frame- 
work but at each iteration, the field is taken as a linear combination of the 
field given by the standard algorithm and its filtered version. The parameter 
of the linear combination is computed numerically to ensure the monotonic 
behavior of the algorithm. 

This algorithm has the advantage of simplicity and general applicability 
whatever the filter chosen. We have presented two examples on molecular 
alignment using either bandpass filters or a discretization in the frequency 
domain mimicking pulse shaping techniques. This work leads to important 
insights into the different ways to achieve molecular alignment. Since there 
exists no unique optimal solution, we can select using spectral filters con- 
trol fields taking into account experimental constraints. For instance, we 
have shown that the control can find a pathway using only three rotational 
frequencies to reach the aligned state. The spectral constraints have sup- 
pressed the other pathways simplifying therefore the structure of the control 
field. 

This algorithm allows to test the choice of the filter with respect to the 
optimal field considered. Indeed, a value of fi^ close to 1 means that the 
filter is not well suited to the form of the electric field determined by the 
monotonic algorithm. This remark is particularly true for the first itera- 
tions of the algorithm since after the first iterations, the optimal field and 
the filtered field can be very close to each other. In the case of pulse shaping 
techniques, we can then determine from the algorithm the minimum number 
of pixels of the mask needed to reach a high alignment efficiency This ren- 
ders the use of optimal control theory more interesting from an experimental 
point of view by making a direct link with control experiments. 
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A Proof of the monotonic behavior of the algo- 
rithm 



In this section, we show the monotonic behavior of the algorithm in the case 
of a nonlinear interaction with the control field. Let E k be the optimal field 
at step k of the algorithm. We determine the field E k +\ at step k + 1 so 
that the variation AJ = J(E k+ \) — J(E k ) > [27]. This variation is given 
by: 

(ioaj = K0/i^ + i(t/))i 2 - \(<f>f\Mtf))\ 2 - r a ^+i(*) 4 - E k {tf\dt. 

Jo 

We introduce the function P k+ i(t) = |{xfc(i)|V'ft+i(*))| 2 - Differentiating it 
with respect to time produces 

(11) Jt Pk+l{t) = ak > k+1 ^t ~ E l+i) 
where 

(12) a Kk+l = 2Im[(^ +1 (t)| X fc(t))(Xfe(i)l^(cos 2 0Aa + a ± )|^+i(t)), 
and a direct integration gives 

(13) P k+1 (t f ) = P k+1 (0) + J*' ^^dt. 

Using the relation |(0/|Vfc+i (*/))| 2 - |(<MV>fc(i/)>| 2 = " flk(*/), 

straightforward computations lead to 

(14) J(£ fc+1 ) - J(E k ) = - I'' [\(t)(E£ +1 - E%) - (El - E 2 k+1 )a k , k+1 ]dt. 

Jo 

The integrand V of J(E k+ i) — J(E k ) can be written as follows: 

(J5)= (E k+1 — E k )[- X(El +1 + El +1 E k + E k+ iEl + El) — a ktk+ i(E k+ i + E k )] . 

We define the optimal electric field E k+ i as the solution of: 

(Mj) + i —E k = rj[—\(E k+1 +El +l E k + E k+ iEl +E k ) — a kjk+ \ (E k+ i +E k )\, 

where rj is a positive constant. It is straightforward to check that J(E k+ i) — 
J{E k ) > for this choice of E k+ \. 
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